En esta sección se presenta una implementación en Python para calcular los parámetros de ecuaciones de estado cúbicas (SRK, PR, RKPR). Las 2 primeras ecuaciónes de estado SRK y PR, son ecuaciones clásicas y ampliamente utilizadas por la industria y la academia que cuentan con 2 parámetros (parámetro de atracción $a_C$ y de repulsión $b$) para describir el comportamiento de sustancias. Por otro lado, la ecuación de estado RKPR, la cual es una propuesta de una ecuación con un tercer parámetro $\delta_1$,el cual permite incluir el efecto estructural de la molécula de la sustancia a la que se quiere describir su comportamiento termodinámico.
Como se mncionó anteriormente en el caso de las ecuaciones de estado (SRK y PR) se tienen los parámetro de atracción $a_C$ y de repulsión $b$ que pueden ser calculados por medio de expresiones que relacionan constantes como temperatura crítica $T_c$, presión crítica $P_c$, volumen crítico $V_c$ y factor acéntrico $\omega$ de una sustancia pura ademas de la constate universal de los gases R.
En el caso de especificar las constantes $T_c$, $P_c$, $V_c$ y $\omega$, es simple calcular los parámetros $a_c$, $b$ y $m$ por medio de las siguientes ecuaciones:
| Parámetros Ecuación de estado SRK | Parámetros Ecuación de estado PR |
|---|---|
| $ a_c = 0.077796070 \frac{R^2 T_c^2} {P_c}$ | $ a_c = 0.45723553 \frac{R^2 T_c^2} {P_c} $ |
| $ b_c = 0.086640 \frac{ R T_c} {P_c}$ | $ b_c = 0.077796070 \frac{R T_c} {P_c} $ |
| $ m = 0.480 + 1.574 \omega - 0.175 \omega^2$ | $m = 0.37464 + 1.54226 \omega - 0.26992 \omega ^2$ |
Ahora en el caso realizar una especificación para los valores de los parámetro de atracción $a_C$, de repulsión $b$ y $m$ para una sustancia pura, es simple obtener los valores correspondientes de las constantes $T_c$, $P_c$, $V_c$ y $\omega$
$$ T_c = \frac{\omega_b a_c} {\omega_a R b} $$$$ P_c = \frac{\omega_b R T_c} {b} $$$$ V_c = \frac{Z_c R T_c} {P_c} $$En el caso del $ \omega$, se debe resolver una ecuación cuadratica que depende de unos determinados valores de constantes $c$, del parámetro $\delta_1$ y el parámetro $m$, que toman determinados valores para casa ecuación de estado:
$$ \omega = 0.5 \frac{- c_2 + \sqrt{c_2^2 - 4 c_1 c_3}}{2c_1} $$| Ecuación de estado SRK | Ecuación de estado PR |
|---|---|
| $\delta_1 = 1.0$ | $\delta_1 = 1.0 + \sqrt{2.0}$ |
| $c_1 = -0.175$ | $c_1 = -0.26992$ |
| $c_2 = 1.574$ | $c_2 = 1.54226$ |
| $c_3 = 0.48 - m$ | $c_3 = 0.37464 - m$ |
En el caso de la ecuación de estado RKPR, se tiene una posibilidad adicional por medio de utilizar el parámetro estructural $\delta_1$ para correlacionar el valor del factor de compresibilidad $Z_c$ con el támaño de la molécula de la sustancia pura que se está tratando. De esta forma, a las especificaciones que se pueden hacer con las ecuaciones de estado (SRK y PR), como en el caso de especificar las constantes como temperatura crítica $T_c$, presión crítica $P_c$, volumen crítico $V_c$ y factor acéntrico $\omega$ de una sustancia pura, se tiene 3 posibles situaciones adicionales:
La primera especificación corresponde a un valor del factor de compresibilidad crítico $Z_c$ y luego determinar el valor del parámetro $\delta_1$ que cumple con dicha especificación. Después se procede con el cálculo del parámetro $k$.
La segunda especificación corresponde a un valor del parámetro $\delta_1$ para el posterior cálculo del parámetro $k$.
La tercera opción es utilizar una correlación termodinámica para obtener un valor de la densidad de líquido saturado de una sustancia pura $\rho(T)_{sat}^{liq}$ y pasarlo como una especificación, para encontrar un valor de los parámetros $\delta_1$ y $k$, que cumplan con la imposición del valor de la densidad de líquido saturado.
Figura 1. Diagrama conceptual del calculo de parámetros ecuación RKPREn la figura 1, los casos de Mode = 1, 2 y 3 corresponden a especificar las constantes ($T_c$, $P_c$, $\omega$) + alguna de las variables ($V_c$, $\delta_1$ , $\rho(T)_{sat}^{liq}$), mientras que el Mode = 4 se refiere a la especificación de los parámetros ($a_c$, $b$, $k$, $\delta_1$) y obtener el valor de las constantes ($T_c$, $P_c$, $V_c$, $\omega$). Este último cálculo se realiza de forma directa como en el caso de las ecuaciones (SRK y PR), por tanto, la siguiente breve explicación se centra en las 3 primeras opciones.
La primera especificación corresponde a dar un valor del parámetro $\delta_1$, con est valor se calcula el factor de compresiblidad $Z_c$ por mediod e las siguientes ecuaciones:
$$d_1 = (1 + \delta_1 ^2) / (1 + \delta_1)$$$$y = 1 + (2 (1 + \delta_1) ^ \frac{1} {3} + \left (\frac{4} {1 + \delta_1} \right)^ \frac{1} {3}$$$$ \omega_a = \frac{(3 y ^2 + 3 y d_1 + d_1 ^ 2 + d_1 - 1)} {(3 y + d_1 - 1) ^ 2} $$$$ \omega_b = \frac{1} {3 y + d_1 - 1} $$$$ Z_c = \frac{y} {3 y + d_1 - 1} $$en $A_0$
factor de compresibilidad crítico $Z_c$, determinado por las constantes ($T_c$, $P_c$, $V_c$):
$$ Z_c = \frac{P_c V_c}{R T_c}$$para el posterior cálculo del parámetro $k$.
La segunda especificación corresponde a un valor del factor de compresibilidad crítico $Z_c$, determinado por las constantes ($T_c$, $P_c$, $V_c$):
$$ Z_c = \frac{P_c V_c}{R T_c}$$que para luego determinar el correspondiente valor del parámetro d1 que cumple con dicha especificación. Después se procede con el cálculo del parámetro k.
La tercera opción es la especificación de un valor para la densidad de líquido saturado a una temperatura. En este caso se debe encontrar un valor d1 y el correspondiente valor del parámetro k, que permita cumplir con la imposición del valor de la densidad de líquido saturado. En este caso, se puede utilizar la clase Thermodynamic_correlations() para obtener un valor para la densidad de líquido saturado de una sustancia pura a una determinada temperatura y luego pasar este valor, como una especificación en la obtención de los parámetros de la ecuación RKPR.
En la figura 4 se muestran como variables de entrada las constantes Tc,Pc, w y alfa que puede ser una especificación de alguno de los 3 parámetros siguientes $(\delta_1, V_c, \rho(T)_{sat}^{liq})$.
La función F1 corresponde a la estimación de un valor para el parámetro {d1} de acuerdo a una correlación preestablecida en el caso de {alfa = Vc}.
La función F2 es el cálculo de los parámetros {ac} y {b} para el correspondiente valor del parámetro {d1}. En el caso de especificar el parámetro {alfa=d1}, el cálculo de los parámetros Zc, ac y b son directos y no requieren de iteración. Mientras que en el caso de alfa = {Vc} se requiere encontrar de forma iterativa el valor del parámetro d1 que verifique el valor de Zc correspondiente por medio del Vc previamente especificado. De manera similar se procede en el caso de {alfa = rho(T)_sat_liq}.
En los ejemplos siguientes se utilizan los datos termodísicos de la base de datos DPPR. Para el caso se tiene como especificación la ecuación de estado RKPR y las constantes criticas para el componente 3-METHYLHEPTANE.
2.3
In [1]:
import numpy as np
import pandas as pd
import pyther as pt
importar las linrerías requeridas, en este caso se trata de las librerías numpy, pandas junto con pyther
In [3]:
properties_data = pt.Data_parse()
component = "3-METHYLHEPTANE"
component = "METHANE"
component = "ETHANE"
component = "PROPANE"
component = "n-HEXATRIACONTANE"
NMODEL = "RKPR"
NMODEL = "PR"
ICALC = "constants_eps"
properties_component = properties_data.selec_component(component)
pt.print_properties_component(component, properties_component)
dinputs = np.array([properties_component[1]['Tc'], properties_component[1]['Pc'],
properties_component[1]['Omega'], properties_component[1]['Vc']])
component_eos = pt.models_eos_cal(NMODEL, ICALC, dinputs)
#ac = component_eos[0]
print(component_eos)
In [6]:
components = ["ISOBUTANE", "CARBON DIOXIDE", 'METHANE', "ETHANE", "3-METHYLHEPTANE", "n-PENTACOSANE",
"NAPHTHALENE", "m-ETHYLTOLUENE", "2-METHYL-1-HEXENE"]
NMODEL = "PR"
ICALC = "constants_eps"
component_eos_list = np.zeros((len(components),4))
for index, component in enumerate(components):
properties_component = properties_data.selec_component(component)
#pt.print_properties_component(component, properties_component)
dinputs = np.array([properties_component[1]['Tc'], properties_component[1]['Pc'],
properties_component[1]['Omega'], properties_component[1]['Vc']])
component_eos = pt.models_eos_cal(NMODEL, ICALC, dinputs)
component_eos_list[index] = component_eos
#components_table_PR = pd.DataFrame(component_eos_list, index=components, columns=['ac', 'b', 'rm', 'del1'])
#print(components_table_PR)
In [7]:
components_table_PR = pd.DataFrame(component_eos_list, index=components, columns=['ac', 'b', 'rm', 'del1'])
print(components_table_PR)
components_table_PR
Out[7]:
In [12]:
name = properties_data.read_dppr()[0]
name
In [ ]:
In [ ]:
De esta forma se observa el calculo simple de los parámetros para la sustancia pura 3-METHYLHEPTANE_RKPR
A continuación se realiza el mismo tipo de calculo pero tomando una serie de 9 sustancias puras, que se pueden extender facilmente a n sustancias, para obtener sus parámetros de nuevo con la ecuación de estado RKPR.
In [3]:
properties_data = pt.Data_parse()
components = ["ISOBUTANE", "CARBON DIOXIDE", 'METHANE', "ETHANE", "3-METHYLHEPTANE", "n-PENTACOSANE",
"NAPHTHALENE", "m-ETHYLTOLUENE", "2-METHYL-1-HEXENE"]
NMODEL = "RKPR"
ICALC = "constants_eps"
component_eos_list = np.zeros((len(components),4))
for index, component in enumerate(components):
properties_component = properties_data.selec_component(component)
pt.print_properties_component(component, properties_component)
dinputs = np.array([properties_component[1]['Tc'], properties_component[1]['Pc'],
properties_component[1]['Omega'], properties_component[1]['Vc']])
component_eos = pt.models_eos_cal(NMODEL, ICALC, dinputs)
component_eos_list[index] = component_eos
components_table = pd.DataFrame(component_eos_list, index=components, columns=['ac', 'b', 'rm', 'del1'])
print(components_table)
Como se observa, los resultados obtenidos son organizados en un DataFrame permitiendo agilizar la manipulación de los datos de una serie de sustancias puras.
In [8]:
components_table
Out[8]:
En el siguiente ejemplo se utiliza la ecuación RKPR pero esta vez con la especificación de la temperatura y densidad de líquido saturado para el CARBON DIOXIDE y de esta forma encontrar el valor del parámetro delta que verifica la especificación realizada para la densidad de líquido saturado.
In [29]:
properties_data = pt.Data_parse()
dppr_file = "PureFull.xls"
component = "CARBON DIOXIDE"
NMODEL = "RKPR"
ICALC = "density"
properties_component = properties_data.selec_component(dppr_file, component)
pt.print_properties_component(component, properties_component)
#dinputs = np.array([properties_component[1]['Tc'], properties_component[1]['Pc'],
# properties_component[1]['Omega'], properties_component[1]['Vc']])
T_especific = 270.0
RHOLSat_esp = 21.4626
# valor initial of delta_1
delta_1 = 1.5
dinputs = np.array([properties_component[1]['Tc'], properties_component[1]['Pc'],
properties_component[1]['Omega'], delta_1, T_especific, RHOLSat_esp])
component_eos = pt.models_eos_cal(NMODEL, ICALC, dinputs)
print(component_eos)